#include <dirent.h>
#include <iostream>
#include <locale.h>
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>


#include "stripe.h"

using namespace cv;
using namespace std;
/*********************************************************************************************************************************************************
*函数描述： stripeException
计算并是否含有横纹的函数，ave求平均值，variance求方差 *函数参数： img
需要计算的图片，BGR存放格式，彩色（3通道），灰度图无效。 *原理：
    1.将图片转成灰度图
    2.求出傅里叶变换后的图
    3.将变换后的图分别投影到x和y轴上
    4.求出x，y轴上的方差，并得出比值。
*函数返回值： 根据比值判断是否存在条纹，无显式返回值
**********************************************************************************************************************************************************/
double ave(double a[], int n) {
  int sum = 0;
  for (int i = 0; i < n; i++)
    sum += a[i];
  return sum / n;
} //平均值，没什么好说的?
  //方差要用的是平均值，不是中值；

double variance(double a[], int n) {
  double sum = 0;
  double average =
      ave(a, n); //函数调用！！不许有【】！！！不许有int和double ！！
  for (int i = 0; i < n; i++)
    sum += (a[i] - average) * (a[i] - average);
  return sum / n;
}
//这个是求方差的，最后要除以n，别忘了！！！

bool stripeException(Mat img) {
  Mat gray, temp;
  cvtColor(img, img, CV_BGR2GRAY);

  int M = getOptimalDFTSize(img.rows); // 获得最佳DFT尺寸，为2的次方
  int N = getOptimalDFTSize(img.cols); //同上
  Mat padded;
  copyMakeBorder(img, padded, 0, M - img.rows, 0, N - img.cols, BORDER_CONSTANT,
                 Scalar::all(0)); // opencv中的边界扩展函数，提供多种方式扩展

  Mat planes[] = {
      Mat_<float>(padded),
      Mat::zeros(padded.size(),
                 CV_32F)}; // Mat 数组，第一个为扩展后的图像，一个为空图像，
  Mat complexImg;
  merge(planes, 2, complexImg); // 合并成一个Mat

  dft(complexImg, complexImg); // FFT变换， dft需要一个2通道的Mat

  // compute log(1 + sqrt(Re(DFT(img))**2 + Im(DFT(img))**2))
  split(complexImg,
        planes); //分离通道， planes[0] 为实数部分，planes[1]为虚数部分
  magnitude(planes[0], planes[1], planes[0]); // 求模，计算幅值
  Mat mag = planes[0];
  mag += Scalar::all(1);
  log(mag, mag); // 模的对数

  // crop the spectrum, if it has an odd number of rows or columns
  mag = mag(Rect(0, 0, mag.cols & -2, mag.rows & -2)); //保证偶数的边长
  // imshow("spectrum magnitude1", mag);
  // waitKey();
  int cx = mag.cols / 2;
  int cy = mag.rows / 2;

  // rearrange the quadrants of Fourier image
  // //对傅立叶变换的图像进行重排，4个区块，从左到右，从上到下 分别为q0, q1, q2,
  // q3 so that the origin is at the image center     //  对调q0和q3, q1和q2
  Mat tmp;
  Mat q0(mag, Rect(0, 0, cx, cy));
  Mat q1(mag, Rect(cx, 0, cx, cy));
  Mat q2(mag, Rect(0, cy, cx, cy));
  Mat q3(mag, Rect(cx, cy, cx, cy));

  q0.copyTo(tmp);
  q3.copyTo(q0);
  tmp.copyTo(q3);

  q1.copyTo(tmp);
  q2.copyTo(q1);
  tmp.copyTo(q2);

  // normalize(mag, mag, 0, 1, CV_MINMAX);  // 规范化值到 0~1 显示图片的需要
  // cout<<mag<<mag.size()<<endl;
  // cout<<"通道数"<<mag.channels()<<endl;

  // imshow("spectrum magnitude", mag);
  // waitKey();
  double W[mag.cols];
  for (int i = 0; i < mag.cols; i++)
    W[i] = 0;
  double H[mag.rows];
  for (int i = 0; i < mag.rows; i++)
    H[i] = 0;
  float count = 0;
  // cout<<"通道数"<<mag.channels()<<endl;
  for (int i = 0; i < mag.rows; i++) {
    for (int j = 0; j < mag.cols; j++) {
      if (mag.at<uchar>(i, j) > 204) {
        // cout<<"每一个点是："<<mag.at<uchar>(i,j)<<endl;
        W[j] = W[j] + 1;
        H[i] = H[i] + 1;
        count++;
      }
    }
  }

  double var1 = variance(W, mag.cols);
  // cout<<"var1是:"<<var1<<endl;
  double var2 = variance(H, mag.rows);
  // cout<<"var2是:"<<var2<<endl;
  double vv = var1 / var2;
  normalize(mag, mag, 0, 1, CV_MINMAX); // 规范化值到 0~1 显示图片的需要
  // imshow("spectrum magnitude", mag);
  // waitKey();
  // cout<<"当前异常亮点数是:"<<count<<endl;
  // count = count/mag.rows/mag.cols;
  if (vv > 15) {
    cout << "条纹密集！！！条纹率是:" << vv << endl;
    return true;
  } else
    cout << "画面正常，条纹率是:" << vv << endl;
  return false;
}
